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Abstract. In this note, we focus on a selection model problem: a mono- 
exponential model versus a bi-exponential one. This is done in the biological 
context of living cells, where small data are available. Classical statistics are 
revisited to improve existing results. Some unavoidable limits are also pointed 
out. 



1. Introduction 

The measurement of molecular dynamic interactions and their respective pro- 
portions in living cells or tissues is a major question in biological and medicine 
research. The Forster resonance energy transfer (FRET) is one of the best known 
approaches to observe and quantitatively study protein-protein interactions at a 
subcellular level ([7]). The FRET measurement can be currently performed by 
fluorescence lifetimes imaging microscopy (FLIM for short) in living cells and tis- 
sus. It can be achieved via the time correlated single photon counting (TCSPC) 
method which provides a lifetime decay curve per site ([B])- To be interpreted, this 
curve is fitted by selecting the "best" (with respect to a given statistical criterion) 
multi-exponential model. Contrary to a mono-exponential model, a bi-exponential 
one witnesses interaction between two proteins. Our aim is to find, pixel per pixel, 
which of these models is accurate. But one difficulty is that the number of observed 
photons per pixel is small for any statistical treatment in order to preserve the liv- 
ing cell and therefore cannot be increased. An attempt to deal with the problem 
can be found in [?]• Our aim here is to go further in this direction pointing out 
some improvements and limits. Some account of statistical methods in this area 
can be found in [5] and [6]. 

1.1. Modelling fluorescence lifetimes. It is not necessary to describe here in 
details FLIM and TCSCP. We only need to understand that lifetimes are measured 
as differences between excitation times (pulses) and emission times of photons. 
Denote by r the period between two consecutive pulses. Here r is 12 nanoseconds, 
near values taken in practice. What is actually measured is a lifetime modulo r 
since we cannot be sure from what pulse it goes. 

It is assumed that lifetimes come from say K species and are observed in the 
interval [0, r) after infinitely many pulses. In these conditions, each lifetime species 
k {1 < k < K) admits the following probability density: 

(1) /fe(t) = afcexp(-a,i)- 



1 ~ exp(-afer) 
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where ak is the inverse mean hfetime of the k-th species. A uniform noise is added 
with density 

(2) /„(t) = irf^. 

r 

If TTfe denotes the proportion of the fc-th species (ttq refers to the noise's one), we 
get the probabihty density of the fluorescence hfetime by writing 

K 

(3) g{t) = Y,^kfkit). 

k=0 

1.2. Modelling the photon emission. Let Ik be the mean photon number of 
species k detected between two pulses. Assume that photons occurrences are inde- 
pendent. Then the total number of detected photons is Poisson distributed with 
intensity TJ2k=o^k if observations take place during T pulses. For a later use, it 
is convenient to set 

(4) ^ = E4. 

fc=o 

Note that we have 

h 

TTk = . 

/ 

Since the noise intensity Iq will be supposed known, it is convenient to consider 
proportions tt^ among all species with fc > 1 except k = 0. Thus, we have for k > I, 

h I 



I-h I -If 



1.3. Maximum likelihood estimation (MLE) and likelihood ratio test. The 

aim is firstly the determination of the most probable parameter 9 := (ai, . . . , a^f , /i, . . 
from observed lifetimes modulo r denoted by ti, . . . ,tn- The noise intensity Iq is 
supposed known. The related log-likelihood is then 

n 

(5) m=mti,.-.,tn)= -/r + nlog(/r)-log(n!)+^log(g(t,)). 

i=l 

For physical reasons, in particular since lifetimes are sure to be between 30pc and 
30ns, we may and do assume that 9 lies in a compact parameter set. 

Numerical optimisation of the likelihood ([5]) is made easier by knowing deriva- 
tives: 



dh 



k 



1=0 

^ = ]l[o,^)(t)^- ^(l-afet-l-Kt-afcr-l)e-"^'-); 

oak ^ (1 - 0-°'='') 

dm = , " , ^MU)-9{U) 

dh 



^ fk{U)-9{U) ^ rj. fk{ti) 



am 



k \ ^ dak 



. . aft,:) / ^ 



9o^k j-{ g{U) I ^ g{ti) 

Denote by 0^ the most probable parameter if there is K species. 

To decide next which model from if = 1 or iiT = 2 is the most accurate, a 
classical statistic is the likelihood ratio 

D:= [m)-m)]- 
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From a theoretical point of view, since we are dealing with the number of compo- 
nents of a mixture model, even the asymptotics under the null hypothesis are not 
the usual statistics. It can be expressed as a supremum over a Gaussian process 
on a subset of a four-dimensional unit sphere (in our case) endowed with the "right" 
covariance function ([T] and references therein). However this process depends also 
on the "true" point 9. Since all calculations are complicated, it is easier to simply 
simulate if we want to know the level of a test associated with a given threshold. 

Notice on the other hand that simulations hint that the likelihood ratio test is 
quite efficient for knowing the number of components in a mixture with compact 
parameter set (see for example U or |2]). 

2. Selection of the number of exponential species K 

2.1. Comparisons. We restricted ourselves to test K = 1 versus K = 2. It can be 
already a difficult and interesting question, if few observed photons are available. 
With the help of simulated observations, we first optimised 9 by MLE for each K 
and next tested K = 1 versus K = 2 via the likelihood ratio statistics D. 

Compared to the one given in [7], the preceding statistical test is as efhcient 
but with about 100 times less observations. For the reader's convenience and for 
comparison, consider the table obtained in [7]: 



Nbr of photons / Ax'^ 


10.0 


20.0 


30.0 


40.0 


50.0 


90.0 


Error(%) 


1000 


35.7 


34.8 


34.3 


34.6 


34.9 


45 


> 20 


10000 


13.7 


12.0 


11.9 


12.1 


12.9 


27.3 


< 20 


100000 


4.2 


1.7 


2.3 


2.7 


4.7 


26.3 


< 5 


1000000 


1.7 


0.0 


0.0 


0.0 


3.3 


32.7 


< 2 



Table 1. Frequency of selection of the wrong model. It depends 
on the observations number and a Ax^ criterion which consists in 
comparing the statistics for K = 1 and K = 2. Simulations 
were performed on a mix of l/ai = 0.6 ns and l/a2 = 2.4 ns with 
different proportions tt[ — 0, .077, .2, .43, 1 with 100 noise photons. 
30 simulations per condition. 



In similar simulation conditions, we have obtained the following: 



Nbr of photons 


Mean error rate (%) 


Best threshold 


Mean error rate 
at threshold 4 


1000 


12.8 


.85 


20 


10000 


0.3 


4 


0.3 


100000 





4 






Table 2. Frequency of selection of the wrong model. It depends 
on the observations number and a likelihood ratio criterion. Simu- 
lations were performed on a mix of 1/ai = 0.6 ns and l/a2 ~ 2.4 
ns with different proportions tt[ = 0, .077, .2, .43, 1 with 100 noise 
photons. 500 simulations for each proportion and number of pho- 
tons. 

Here, mean error rate is the average over the simulation number of the percentage 
to select the wrong model. Best threshold means threshold that gives the smallest 
mean error rate; using a very crude optimisation. Notice that the strange values 
0, .077, .2, .43, 1 of tt'i's proportion correspond to values .25, .5, .75, 1 of proportion 
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r]i = — — — — considered in 7 . A consequence is that we never test more 

short-life photons than long-hfe. Moreover the case tt[ = .077 is not very far away 
from the mono-exponential case. In particular, with 1000 photons, among which 
100 noise photons, the expected number of photons with 0.6 nanoseconds lifetime 
is less than the number of noise photons. If we compute the error rate for 1000 
photons without that case, we obtain for instance a 2.6% error rate for the likelihood 
ratio test at threshold 3. 

2.2. Simulation scheme. The data set generation algorithm is as follows: 

1. Sample rife the number of photons for each species k, including noise (fc — 0), 
from a Poisson distribution of parameter Tlk- 

2. Draw rik lifetimes with distribution density fk for each species k. 

3. Return the set of all the sampled lifetimes, regardless of k. 
Some differences between simulation methods should be noted: 

• We use random Poissonian number of photons rather than fixed number of 
photons : we take into account the "offset noise" . 

• Instrumental response: we neglect the .03 nanoseconds long instrumental 
response function. 

• Exact times vs channels: we did not use bins and worked as if we knew the 
exact detection times. 

Nevertheless these differences should have little effect and comparisons still make 
sense. 

3. Further comments 

3.1. With closer lifetimes. If we choose 1/ai = 1 ns and l/a2 = 2 ns as mean 
lifetimes, it is harder to select the right number of species: 

• With 10000 photons and — = .01 as noise ratio: 

1 - TTq 

— If TT^^ = = .5 or tt[ — .75, ttj — .25, no wrong selection should occur, 

— li Tr[ = .25, TTj = .75, the error rate is about .1% when the threshold is 
calibrated so as to balance errors "mono towards bi" and "bi towards 
mono" . 

• With 1000 photons and = .01 as noise ratio: if t^'i = '^'2 = -5, the 

1 - TTo 

error rate is about 15% when the threshold is calibrated so as to balance 

errors "mono towards bi" and "bi towards mono" . 
If we choose close mean lifetimes such as 1/ai — 1.4 ns and l/a2 = 1-6 ns, we are 
too close to the "border" of the model, and about 1 million photons is required 
to distinguish the two components. By border, we mean a proportion close to 
or 1/ai close to l/a2 so that identifiability problems occur with small samples. 
Asymptotically, when we get n times closer to the border of a mixture model, 
we need times as many photons to get the same statistical efficiency, for any 
procedure [2]. 

3.2. Absolute limits. The former sentence about rates when we get nearer the 
border is a first expression of limits that cannot be broken, no matter how smart 
the statistical procedure. To give a small taste of what to expect, here are the best 
error rate when having to choose specifically between two possible sets of lifetime 
parameters and corresponding distribution probabilities and with equal a 
priori probabilities. In that situation, which is easier than the one studied in the 
article, the optimal choice is the one with greater observed likelihood, and the error 
rate is 1 - i ||/i - 
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• With 32 observed photons and a signal to noise ratio of 1/10, choose be- 
tween a mono-exponential with lifetime 2.4 ns, and a bi-exponential with 
proportions 0.077 and 0.923 and lifetimes 0.6 and 2.4 ns: optimal error rate 
> 25%. 

• With 32 observed photons and no noise, choose between a mono-exponential 
with lifetime 2.6 ns, and a bi-exponential with proportions one half and 
lifetimes 2.5 and 2.7 ns: optimal error rate > 49.75%. 

The second case is almost as bad as a coin toss, ignoring the data. 
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